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A  numerical  study  on  laminar  incompressible  flows  in  2-D  straight  walled 
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downstream  the  diffusers.  Geometries  with  area  ratios  AR=1.15  to  5  and  non- 
dimensional  lengths  of  L/Wi=l  to  48  are  considered.  Results  are  presented  in  terms  of 
flow  regime  maps  for  Reynolds  numbers  of  105,  210,  314,  420,  629,  1,048  and  pressure 
recovery  coefficients  maps  for  Re  numbers  of  105,  210,  314,  420  and  629.  In  addition 
time  resolved  simulations  of  impulsively  starting  flow  are  considered  at  Re=210,  314  for 
12  geometries  on  the  flow  regime  map.  Four  flow  regimes  can  be  distinguished 
depending  on  diffuser  geometry.  With  increasing  divergence  angle  the  flow  goes  from 
attached  to  symmetrically  separated  to  asymmetrically  separated  and  finally  to  a  non  2-D 
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I.  INTRODUCTION 


A.  OVERVIEW 

A  diffuser  is  a  common  device  used  in  internal  flow  systems.  Such  flow  systems 
are  encountered  in  turbomachines  between  a  compressor  and  a  combustor  or  at  the  exit  of 
a  turbine,  duct  flows,  flow  meters,  aircraft  engine  inlets,  closed  circuit  wind  tunnels,  and 
pumps.  Diffusing  passages  may  also  be  encountered  in  small-scale  devices  such  as  fluidic 
actuators  and  Micro-Electro-Mechanical  Systems  (MEMS). 

The  primary  purpose  of  using  a  diffuser  in  most  of  the  aforementioned 
applications  is  to  maximize  static  pressure  recovery  while  minimizing  total  pressure  loss 
along  the  direction  of  the  flow.  In  other  words  a  diffuser  serves  as  a  converter  of  the  flow 
dynamic  head  to  static  pressure.  This  is  achieved  by  small  angle  divergent  walls,  in  which 
the  flow  is  confined,  resulting  in  a  certain  cross  sectional  area  increase. 

Unlike  nozzles  where  the  static  pressure  decreases  in  the  flow  direction  (favorable 
pressure  gradient),  the  fluid  particle  flowing  through  a  diffuser  experiences  an  increasing 
static  pressure  (adverse  pressure  gradient)  resulting  to  a  number  of  fluid  dynamic 
phenomena.  These  phenomena  range  from  flow  separation  and  unsteadiness  to  transitory 
stall  and  violent  flow-excited  static  pressure  fluctuations.  Most  of  these  flow  situations 
are  very  complicated  and  although  they  are  qualitatively  understood  they  are  hard  to 
predict  quantitatively. 

In  Figure  1,  2,  and  3  the  most  common  and  simple  diffuser  shapes  and  relevant 
geometric  parameters  are  shown.  However  one  can  have  more  general  shapes  of  diffusers 
like  the  curved  ones  shown  in  Fig.4,  which  depicts  diffusing  passages  in  an  axial 
compressor  stage. 
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Figure  3.  Annular  diffuser 


Figure  4.  Stator  blades  diffuser. 

Because  of  the  complexity  of  the  diffuser  problem  a  great  number  of  different 
variables  must  be  considered  in  order  to  study  diffuser  flows.  The  parameters  driving 
diffuser  flow  behavior  and  performance  can  be  divided  in  two  categories.  The  first  is 
diffuser  geometrical  shape  parameters.  The  second  is  flow  parameters. 
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The  geometric  parameters  are  generally  described  as  non-dimensional  length 
L/Wi  which  is  the  ratio  of  the  diffuser  length  to  throat  width,  the  area  ratio  between  outlet 
and  inlet  cross  sections  AR,  the  half  angle  of  divergence  0,  the  inlet  ratio  of  the  radii  in 
the  case  of  annular  diffusers  or  the  aspect  ratio  AS  (width  to  height)  in  the  case  of 
rectangular  diffusers,  depending  on  the  choice  made  and  the  geometry  in  hand. 

The  flow  parameters  are  related  to  the  flow  conditions  in  hand.  Inlet  flow 
variables  such  as  turbulence  intensity,  inlet  swirl,  boundary  layer  thickness  (blockage 
factor),  shape  of  the  velocity  profile  at  the  inlet,  presence  of  wakes  or  stall  regions  at  the 
inlet,  throat  Reynolds  number,  throat  Mach  number,  inlet-outlet  flow  blockage,  other 
flow  devices  situated  upstream  or  downstream  of  the  diffuser,  are  some  of  the  flow 
parameters  driving  diffuser  performance  and  behavior. 

Typically  three  factors  are  used  to  quantify  efficiency  and  performance  of  a 
diffuser.  The  first  is  a  pressure  recovery  coefficient  Cp ,  defined  as  the  ratio  of  the 

difference  between  exit  static  pressure  pout  and  inlet  static  pressure  pin  to  the  throat 
dynamic  head  and  used  mostly  for  subsonic  flows: 


’  _  Pom  Pin 

p  \llpU 2 


Where  p  is  the  density  and  U  is  the  mass  averaged  flow  velocity.  The  second 
factor  is  a  loss  coefficient  nd  defined  as  the  ratio  of  the  outlet  to  inlet  total  pressures 
PTout ,  PTin  and  used  for  supersonic  flows  where  the  total  pressure  loss  in  the  diffuser  is 
more  significant: 


= 


The  third  factor  is  diffuser  effectiveness  s .  It  is  defined  as  the  ratio  of  measured 
pressure  recovery  C  to  ideal  inviscid  pressure  recovery  C  . : 
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In  general  it  is  desired  to  optimize  diffuser  design  by  maximizing  the  above 
factors.  High  Cp  means  that  the  diffuser  recovers  the  greatest  amount  of  static  pressure 

possible  from  the  incoming  dynamic  head.  High  s  means  a  good  static  pressure  recovery 
relative  to  a  maximum  possible  ideal  value  for  inviscid  flow  through  the  same  diffuser. 
However  there  are  applications  where  other  characteristics  of  the  flow  field  are  of  interest 
apart  from  max  recovery  or  effectiveness.  It  may  be  needed  to  have  as  uniform  a  flow  as 
possible  at  the  outlet  of  the  diffuser  or  flow  free  of  pressure  oscillations  (i.e.  combustion 
flows). 

In  many  cases  diffuser  performance  augmentation  takes  place  by  some  means  of 
separation  control  such  as  wall  suction  or  blowing,  vortex  generators,  vanes,  screens,  etc. 
The  objective  is  to  reenergize  the  momentum  deficient  flow  near  the  diffuser  walls  such 
that  separation  does  not  take  place. 

B.  LITERATURE  REVIEW  ON  SUBSONIC  DIFFUSERS 

Research  work  on  diffuser  flow  spans  several  decades.  Gibson  G.H.  and  Eiffel  G. 
reported  experimental  results  as  early  as  1910  and  1919  respectively  [23,24],  More 
attention  was  drawn  to  diffuser  flows  in  the  early  fifties  at  the  same  time  with  the  advent 
of  high-speed  flight  and  the  increase  in  need  to  have  more  efficient  and  reliable  internal 
combustion  engines  in  which  the  diffuser  is  an  important  integral  part.  As  a  result  high 
Reynolds  number  subsonic  and  supersonic  flows  were  considered  in  the  early 
experimental  investigations.  Because  this  work  will  consider  only  subsonic  flows,  a  brief 
review  of  prior  work  in  such  flows  will  be  presented  here.  The  various  parameters  and 
their  influence  on  diffuser  performance  will  be  explained.  Such  prior  work  on  subsonic 
diffusers  can  be  divided  as  follows: 

1)  Experimental/analytical  work  on  steady  incompressible  flow  in  diffusers  of 
various  geometric  shapes.  The  flow  at  the  diffuser  throat  has  been  taken  to  have  an 
inviscid  core  with  thin  turbulent  or  laminar  boundary  layers. 

2)  Experimental/analytical  work  on  unsteady  flow  though  diffusers  with  throat 

flow  characteristics  same  as  the  ones  mentioned  above  in  1 . 
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3)  Experimental  work  on  improving  performance  by  controlling  the  diffusing 
flow  separation  using  various  means  (vortex  generators,  wall  jets,  wall  suction,  screens  or 
vane  insertion,  and  others). 

4)  Numerical  methods  applied  to  compute  diffuser  flows. 

In  1959,  Kline  et  al  [1,  4]  reported  data  on  performance  and  flow  regimes  for 
incompressible  high  Reynolds  number  flow  (Re>70,000)  in  two  dimensional  straight  wall 
diffusers.  The  authors  provided  and  an  extensive  literature  survey.  In  their  experiment 
they  varied  diffuser  length  to  inlet  height  L/Wi  from  1.5  to  25  and  the  area  ratios  AR 
from  1.3  to  10.  The  inlet  velocity  profiles  used  had  an  inviscid  core  and  small  turbulent 
boundary  layers  close  to  the  walls.  The  data  are  correlated  in  flow  regime  maps  and  the 
following  important  observations  were  made:  The  flow  in  a  diffuser  when  the  throat 
width,  wall  length  and  inlet  flow  conditions  are  held  constant,  while  increasing  the  angle 
of  divergence  from  zero  degrees,  can  be  classified  in  four  distinct  regimes: 

1)  A  region  of  no  appreciable  stall  in  which  the  main  flow  is  almost  steady  and 
not  separated, 

2)  A  region  of  large  transitory  stall  in  which  the  flow  is  unsteady  and  highly 
pulsating, 

3)  A  region  of  fully  developed  stall  on  one  wall  of  the  diffuser  where  the  flow  is 
steady,  and 

4)  A  region  in  which  jet  flow  emerges  from  the  throat  and  separation  begins  just 
down  stream  of  the  throat. 

All  the  above  flow  regimes  were  presented  in  a  map  shown  below  in  Fig.  5.  They 
also  discussed  if  and  which  of  the  flow  parameters  might  have  the  most  influence  on  the 
above  flow  regimes.  One  major  characteristic  is  the  independence  of  the  flow  regimes  on 
Reynolds  number  in  the  range  tested  above.  It  is  also  found  that  increasing  the  free 
stream  turbulence  intensity  has  the  effect  of  shifting  a-a  line  of  no  appreciable  stall  down. 
In  regard  to  the  effect  of  aspect  ratio  AS  on  flow  regimes  the  data  suggested  that  for 
aspect  ratios  above  1 :4  there  is  no  appreciable  change  in  the  flow  regimes. 
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In  a  1961  paper  Waitman  et  al  ([2])  reported  an  experimental  study  on  the 
influence  of  inlet  boundary  layer  thickness  in  the  performance  of  two-dimensional 
subsonic  diffusers.  Several  diffusers  were  tested  spanning  inlet  lengths  to  heights  from 
8.0  to  48.0,  total  divergence  angles  from  2.5  to  40  degrees.  Inlet  boundary  layers  were 
turbulent  and  varied  from  extremely  thin  to  fully  developed  turbulent  channel  flow.  Low 
Mach  number  (<0.2),  high  throat  Reynolds  number  (Re>  10,000)  were  used.  Flow 
obstructions  upstream  of  the  diffuser  throat  were  also  used  for  brief  investigation  on  the 
effects  of  their  wakes  on  performance.  It  was  found  that  pressure  recovery  is  a  strong 
function  of  inlet  boundary  layer  thickness,  decreasing  with  increased  thickness. 
Turbulence  intensity  was  found  to  increase  pressure  recovery  as  intensity  increased.  The 
flow  regimes  chart  given  in  ref.  [1]  was  altered  by  small  degree  in  the  sense  that  stalls 
appeared  in  smaller  value  of  angle  of  divergence  as  boundary  layer  thickened  at  the  inlet. 

In  1965  [3]  Sorvan  and  Klomp  gave  an  extensive  experimental  investigation  of 

pressure  recovery  on  a  large  number  of  annular,  low  subsonic  diffusers  (Ma<0.3).  Also 
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data  for  rectangular  and  conical  diffusers  from  previous  studies,  like  [1],  [2],  were 
gathered  and  used  for  evaluation  of  performance  characteristics  for  all  three  types  of 
geometries.  The  study  showed  a  significant  degree  of  similarity  on  performance  among 
the  three  types  of  diffusers,  if  N IWX,  N  /  Rv  L/  ARt  are  taken  as  the  non-dimensional 

length  parameter  for  rectangular,  conical,  annular  units  respectively.  The  concept  of  area 
blockage  was  shown  to  correlate  the  effect  of  inlet  boundary  layer  thickness  on  pressure 
recovery  for  all  three  types  of  diffusers.  Area  blockage  occurs  in  internal  flow  sytems 
because  of  the  non-uniformities  in  velocity  over  the  cross  section  of  the  systems.  The 
idea  of  “inefficient”  versus  “insufficient”  diffusion  is  introduced;  the  first  describing 
viscous  effects  on  diffuser  effectiveness  while  the  second  implying  deviations  from  ideal 
one-dimensional  inviscid  flow  due  to  velocity  profile  non-uniformities,  ft  is  shown  that 
the  reduction  in  pressure  recovery  due  to  boundary  layer  thickening  is  primarily  the  result 
of  “insufficient”  rather  than  “inefficient”  diffusion. 

Wolf  and  Johnston  [5]  investigated  the  effect  of  non-uniform  inlet  velocity 
profiles  experimentally.  A  two-dimensional  rectangular  diffuser  is  chosen  and  the  inlet 
velocity  profile  is  changed  to  include  wakes,  non-uniform  shear  and  jet  type  flows.  The 
most  important  of  their  results  is  the  increase  in  pressure  recovery  for  several  geometries 
when  a  wake  is  present  in  the  inlet  profile.  However  all  other  types  of  inlet  non¬ 
uniformities  decrease  recovery.  Wake  non-uniformities  are  amplified  downstream  the 
diffuser.  Also  diffuser  flow  regimes  are  shown  to  depend  on  the  relative  location  of  the 
low  velocity  fluid  at  the  inlet  profile.  If  the  low  velocity  fluid  is  near  to  one  inlet  wall, 
separation  tends  to  appear  on  that  wall. 

McMillan  and  Johnston  [6],  [7]  treat  a  low  aspect  ratio  (AS=0.1)  two-dimensional 

diffuser  with  fully  developed  turbulent  inlet  flow.  The  importance  of  AS  alone  is 

investigated  experimentally  in  the  Reynolds  number  range  from  20,000  to  70,000.  The 

specific  feature  of  this  type  of  diffuser  is  the  close  proximity  of  the  parallel  walls  and  the 

effect  that  this  might  have  on  diffusion,  since  the  “presence”  of  these  walls  through 

viscosity  can  be  felt  in  the  core  flow.  The  results  show  that  flow  regimes  are  more 

dependent  upon  inlet  conditions  compared  to  high  aspect  ratio  diffusers.  Regime 

transitions  depend  strongly  on  inlet  Re  number  as  well  as  the  diffuser  geometry.  On  the 

contrary  experience  with  high  aspect  ratio  diffusers  shows  that  flow  regime  depends 
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primarily  on  diffuser  geometry  and  only  slightly  on  inlet  flow  conditions  and  Re  number. 
Pressure  recovery  was  found  to  increase  with  increasing  Reynolds  number  for  the  low 
aspect  ratio  diffuser.  Losses  of  total  pressure  along  the  diffuser  are  comparable  to  those 
of  a  straight  channel  of  the  same  length  as  the  diffuser. 

Runstadler  and  Dolan  [8,  9]  investigate  high  subsonic  flow  in  two-dimensional 
diffusers.  They  employ  geometries  of  different  lengths,  aspect  ratios  from  0.25  to  5,  inlet 
Ma  numbers  from  0.2  to  1.0;  blockage  factors  from  0.02  to  0.12.From  the  collected 
experimental  data  performance  maps  are  constructed  while  aspect  ratio,  Ma  number, 
blockage  and  Re  number  being  held  constant.  The  diffuser  performance  maps  for 
compressible  flow  are  significantly  different  from  the  maps  taken  for  incompressible 
flows.  Higher  blockage  reduces  the  pressure  recovery.  Inlet  Ma  was  shown  to  have  a  very 
mild  increasing  effect  on  recovery  up  until  chocking,  beyond  witch  the  recovery 
deteriorates  rapidly.  Increasing  Reynolds  number  (varied  from  170,000  to  580,000) 
produces  an  increase  in  pressure  recovery. 

Dighe  and  McDonald  [10]  investigated  the  effect  of  aspect  ratio  on  plane-wall 
diffuser  performance  and  flow  regimes  for  Re  numbers  in  the  range  3,000<Re<28,000 
and  with  laminar  boundary  layers  at  the  diffuser  throat.  They  employ  twelve  diffuser 
geometries  with  aspect  ratios  0.21<AS<2.65  and.  All  diffusers  used  had  constant  length 
N/Wi=18  and  total  divergence  4°  half  angle.  The  results  show  that  Re  number  has  a 
profound  effect  on  pressure  recovery  for  all  aspect  ratios  tested,  increasing  Re  number 
from  3,000  to  20,000  increases  recovery  from  0.4  to  0.75  at  constant  aspect  ratio  .  Even 
though  diffuser  geometries  are  chosen  to  be  in  the  unstalled  flow  regime,  when  the  throat 
Re  number  was  low  3000<Re<5000,  flow  regime  changed:  This  is  an  indication  that  for 
lower  Reynolds  number  flows  Reynolds  number  is  a  controlling  flow  regime  parameter. 
Jet  flow  was  observed  for  low  aspect  ratio  diffusers  (AS<0.94),  while  two  dimensional 
stall  was  exhibited  in  diffusers  with  AS>0.94.  Above  Re=5000  none  of  the  diffusers 
investigated  in  this  study  exhibited  stall. 

Stenning  [11],  Smith  [12],  and  Kwong  [13]  investigated  unsteady  diffusers  flows 
in  a  Re  number  range  of  3,000  to  20,000.  They  conclude  that  amplification  of  inlet  mean 
flow  oscillations  in  the  diffuser  take  place  when  throat  Re  number  is  in  the  region  of  3000 
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to  6000,  self-excited  low  frequency  flow  oscillations  occur  in  the  transitory  stall  regime 
(for  Re  numbers  higher  than  12,000-15,000)  and  that  an  amplification  of  such  oscillations 
happens  when  inlet  disturbances  are  present  with  periods  0.5  Tm  to  Tm,  were  Tm  is  the 
mean  oscillation  period  of  the  stall  regions  from  one  wall  of  the  diffuser  to  the  other. 

Various  methods  can  be  used  to  enhance  performance  and  /or  outlet  flow  stability 
and  uniformity.  The  simplest  ways  employ  devices  like  vortex  generators,  vanes  or 
screens,  calling  them  passive  separation  control.  Wall  suction  to  take  away  the  low 
momentum  fluid  near  the  wall  or  blowing  to  add  momentum  can  be  employed.  These 
methods  are  called  active  separation  control.  The  general  idea  used  in  all  methods  is  to 
reenergize  the  boundary  layer  that  forms  on  the  diffuser  walls,  thus  preventing  low 
momentum  fluid  from  separating  due  to  the  adverse  pressure  gradient.  It  can  not  be  said 
that  one  method  is  more  successful  over  the  other  because  depending  on  requirements  of 
the  design  one  might  exploit  advantages  of  a  specific  method  that  are  not  present  in 
others.  See  references  [14],  [15],  [16]  on  diffuser  performance  augmentation  by  flow 
separation  control. 

In  the  recent  years  some  effort  has  been  done  to  compute  diffuser  flow  fields. 
References  [17],  [18],  [19],  [20]  are  papers  devoted  in  investigating  diffuser  flows 
computationally.  Successful  computations  at  high  Re  numbers  (where  modeling  of 
turbulence  is  needed)  are  feasible  if  the  flow  inside  the  diffuser  is  not  separated.  That 
means  very  low  diffusion  angles.  For  separated  flows  only  laminar  very  low  Re  (<1 14) 
number  computations  are  shown  to  be  feasible  [18]. 

In  summary  most  of  the  work  reported  on  subsonic  diffuser  flows  is  experimental 
for  medium  to  high  Reynolds  numbers  (>2,000).  This  work  has  identified  the  parameters 
and  the  extent  each  parameter  affects  diffuser  performance  and  flow  regimes.  It  was 
shown  that  geometric  variables  like  AR,  AS,  L/Wi  and  flow  related  variables  like  Re, 
Ma,  throat  boundary  layer  thickness  and  turbulence  intensity  affect  diffuser  behavior  in 
various  degrees  mentioned  in  the  preceding  paragraphs.  The  results  presented  in  terms  of 
flow  regime  and  pressure  recovery  coefficient  maps  are  used  to  design  optimum  and 
efficient  diffusers  operating  in  the  high  Reynolds  number  range.  However  little  work  is 
reported  for  diffuser  flows  in  the  Re  number  range  below  2,000.  In  this  Reynolds  number 
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range  the  optimum  characteristics  of  the  diffusers  have  not  been  extensively  investigated. 
In  addition  the  effect  of  geometric  and  flow  parameters  on  diffuser  behavior  is  not 
extensively  studied. 

C.  OBJECTIVES  AND  SCOPE  OF  PRESENT  WORK 

The  objectives  of  this  study  were: 

1)  Validate  the  commercial  computational  fluid  dynamics  (CFD)  code  CFD- 
ACE+  against  a  known  analytical  solution  of  the  2-D  Navier-Stokes  equations  for 
straight  plane  diffuser  flows. 

2)  Investigate  different  types  of  flow  regimes  that  might  be  encountered  in  a 
straight  plane  diffuser  with  a  parallel  channel  at  both  inlet  and  exit  while  varying  the  inlet 
Re  number  in  the  range  105-1048  and  the  geometry  of  the  diffuser  in  the  range  of  1.15-5 
for  the  AR  and  1-48  for  the  L/Wi.  For  high  Re  flows,  increasing  diffusion  angle  (all  other 
geometric  parameters  held  constant),  the  flow  undergoes  various  configurations  called 
flow  regimes  [1],  In  this  work  a  parallel  attempt  is  made  but  for  low  Re  number  fully 
developed  flows  and  by  purely  computational  means. 

3)  Determine  optimum  performance  characteristics  in  terms  of  maximum  static 
pressure  recovery  coefficient  for  the  same  diffuser  geometries  and  Re  number  range 
mentioned  in  2  above. 

D.  ORGANIZATION 

Chapter  II  presents  a  brief  description  of  numerical  methods  and  features  of  the 
software  package  used  in  this  work. 

Chapter  III  discuses  exact  solutions  of  the  Navier-Stokes  equations  for  Jeffery- 
Hamel  flows  in  a  2-D  diffuser  and  the  accuracy  of  the  software  against  these  solutions. 

Chapter  IV  describes  the  numerical  solution  of  the  Navier-Stokes  equations  for 
several  diffuser-outlet  channel  configurations  and  Re  numbers  in  the  50-1048  interval 
using  the  commercial  CFD-ACE+  code.  The  results  are  correlated  in  flow  regimes  and 
pressure  coefficient  performance  maps. 
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Chapter  V  presents  the  time  resolved  numerical  calculations  done  in  CFD-ACE+ 
for  Re=320,  210  and  several  diffuser  geometries  in  the  two  most  important  flow  regimes 
found  in  chapter  IV. 

Chapter  VI  concludes. 
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II.  OVERVIEW  OF  CFD-ACE+  AND  THE  SPECIFIC  FEATURES 
OF  THE  PROGRAM  USED  IN  THIS  STUDY 

A.  GENERAL  DESCRIPTION 

The  commercial  program  CFD-ACE+  is  a  solver  that  utilizes  a  finite  volume 
pressure  correction  based  method  to  solve  the  conservation  equations  of  mass, 
momentum  and  energy.  It  provides  the  user  three  graphical  user  interfaces  (GUI’s) 
namely  a  preprocessor  CFD-GEOM  for  geometry  creation,  a  processor  CFD-ACEU 
which  is  the  actual  solver  of  the  descretized  equations  and  a  post  processor  CFD-VIEW 
for  viewing  and  processing  the  results. 

B.  NUMERICAL  METHODS 

1.  Discretization  Schemes 

The  software  provides  several  different  schemes  for  discretizing  the  partial 
differential  equations  in  space  and  time.  Throughout  this  study  the  second  order  scheme 
will  be  used  for  spatial  discretization  and  the  Crank-Nicholson  scheme  for  time  derivative 
discretization.  The  resulting  finite  difference  equations  (FDE)  are  nonlinear  algebraic 
equations  formed  for  each  computational  shell  in  the  following  form: 

(CCp  —  Sp  )<f>p  =  ^  a„b<t)nb  + 

nb 

Where,  (/>  is  the  dependent  variable  and  a  the  non-linear  coefficients.  An 
iterative  method  is  employed  to  solve  the  set  of  the  resulting  system  of  non-linear 
algebraic  equations.  To  formulate  an  equation  for  the  unknown  variable  pressure  the 
continuity  is  used  by  adopting  the  SIMPLEC  algorithm.  For  further  information  on  the 
code  the  reader  is  referred  to  the  CFD-ACEU  user  manual  and  the  additional  references 
therein  [22],  Convergence  is  achieved  if  the  residuals  of  the  dependent  variables 
(pressure,  u  velocity  and  v  velocity  in  this  study)  drop  several  orders  of  magnitude  (>5) 
and  the  continuity  equation  is  satisfied  by  at  least  six  orders  of  magnitude  after  a  number 
of  iterations  has  been  performed.  Note  that  from  this  point  on  the  residual  drop  and  the 
continuity  equation  satisfaction  with  the  number  of  iterations,  will  be  referred  to  as 
“residual  history”  and  “mass  flow  summary”  respectively. 
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2.  Boundary  Conditions 

A  number  of  options  on  boundary  conditions  are  available  to  the  user.  Indicated 
here  are  the  boundary  conditions  that  will  be  used  in  this  study: 

1)  Specification  of  u  velocity  profiles  or  uniform  total  pressure  at  the 
computational  domain  inlet. 

2)  Static  pressure  at  the  computational  domain  outlet. 

3)  No  slip  no  penetration  conditions  at  the  walls. 

3.  Grid  Generation 

The  software  allows  structured  and  unstructured  grid  generation.  In  this  work  the 
structured  grid  is  used.  For  further  information  about  the  meaning  of  each  type  of  grid 
the  reader  is  referred  to  the  manual.  Throughout  this  study  computational  shell  aspect 
ratio  and  skew  ness  are  kept  to  a  minimum.  The  maximum  shell  aspect  ratio  is  kept 
under  8  and  the  smallest  angle  in  any  one  shell  is  kept  under  70°. 
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Figure  6.  Computational  shell  characteristics. 


C.  SPECIFIC  CFD-ACEU  SETTINGS  USED  IN  THIS  STUDY 

In  this  thesis  only  the  “flow”  module  of  the  solver  will  be  used  since  the  flow  is 
assumed  isothermal  and  the  energy  equation  need  not  be  solved.  Inlet  parabolic  velocity 
boundary  conditions  will  be  imported  as  data  files  after  they  are  generated  in  MATLAB. 
Relaxation  for  the  pressure  correction  equation  and  velocity  variables  of  0.1 -0.2  and 
0.05-0.09  respectively  will  be  used. 
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III.  ANALYTICAL  AND  NUMERICAL  SOLUTION  OF  THE 
NAVIER-STOKES  EQUATIONS  IN  A  2-D  DIFFUSER  FOR 
JEFFERY-HAMEL  FLOWS 


A.  PURPOSE 

The  Jeffery-Hamel  hydrodynamic  flow  problem  in  a  plane  diffuser  is  a  well- 
known  exact  analytical  solution  of  the  two-dimensional  steady  incompressible  Navier- 
Stokes  equations  [21].  In  this  study  the  accuracy  of  the  CFD  software  (used  later  in  the 
study)  will  be  checked  against  the  analytical  solution  and  it  will  be  shown  that  the 
multiplicity  of  solutions  predicted  analytically  is  also  observed  in  the  direct  numerical 
solution  of  the  Navier-Stokes  equations  for  flows  computed  in  this  study. 

B.  ANALYTICAL  FORMULATION  AND  SOLUTION  OF  THE  PROBLEM 

1.  Formulation  of  the  Problem 

The  hydrodynamic  problem  is  described  from  the  continuity  and  the  Navier- 
Stokes  equations  written  in  polar  coordinates  form  as: 


d(ru )  dv 


(1) 


dr  d6 


(2) 


dv  v  dv  uv  1  dp  2  du  v . 

u  —  + - +  —  = - — +  v(V  v  +— - -) 

dr  r  dO  r  pr  d6  r2  d6  r2 


dv  v  dv  uv 


r2  dO  r 


The  coordinate  system  is  given  below: 


Figure  7.  Diagram  of  the  2-D  diverging  channel  with  the  line  source  at  the  origin. 
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If  the  flow  is  assumed  to  be  purely  radial  emerging  from  a  line  source  that  is  if 
u  =u(r,0),v  =  0  from  continuity  (1), 


d(m) 


=  0 


ru  =  fnc(0 ) 


u  = 


F{9) 


dr  r 

Substituting  u  =  u(r, 0),v  =  0  back  in  the  momentum  equations  (2)  we  get, 


(3) 


du 

u  —  = 
dr 


1  dp 
— —  +  v 
p  dr 


d2u  1  du  u  1  d2u 
dr2  r  dr  r2  r2  dO2 


n  1  dp  2v  du 

0  = - — +  — r - 

pr  dO  r2  dO 


From  (3)  calculate  all  derivatives  in  (4)  with  respect  to  F: 

du__F_  du___F^  dhi__2F_  d2u  _ 
dr  r2  ’  dO  r  dr2  r 3  ’  dO2  r 

Substitute  in  (4)  to  get 


dp  _  pF 2  pF" 


+  - 


dr  r"  r" 
dp  _  2 pF' 
l  d0~  r3 


Eliminate  pressure 

dp  _  pF2  pF" 


-  +  - 


dr 

dp  _  2 pF' 
dO  r3 


>pF'"  +  2pFF'  +  4pF'  =  0~- 


(4) 


(5) 


(6) 


(7) 


F"'  +  -FF'  +  4F'  =  0 

V 


(8) 


Where:  F  is  the  similarity  velocity  profile.  The  boundary  conditions  in  terms  of 
the  new  function  F  used  to  solve  the  third  order  ode  become: 


F(±a)=0 

(9) 

F\ 0)  =  0 

(10) 
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If  equation  (3)  is  normalized  by  the  value  of  the  velocity  profile  at  0=0  ,  F( 0) 


and  6  by  the  max  half  angle,  a  of  any  given  diffuser,  they  take  the  following  form: 

u(r,0)  =  um 


C(ff\-FW  n-° 

^(O)  a 


FQ 0)-^j+  -  F(0)G  —  +  4F(0)—  =  0 
a  dij  v  adt)  adrj 


Or 


2F(0)a  dG  t  A_2dG 


drf 


-  +  ■ 


v 


G - +  4aA - =  0 

dr)  dr) 


G'"  +  2a  Re  GG'  +  4a  G' = 0 


(11) 

(12) 

(13) 

(14) 

(15) 


Where  G  is  the  normalized  similarity  profile  and  Re  is  a  Reynolds  number 


Re= 


ura 


v 


(16) 


Note  for  small  channel  angles  a,  ra  is  approximately  equal  to  the  channel  half 
height.  The  boundary  conditions  are  transformed  to: 


G(±1)=0 

(17) 

G"(0)=0 

(18) 

G(0)  =  1 

(19) 

This  non-linear  third  order  ode  along  with  the  boundary  conditions  is  solved  for 
the  invariant  velocity  profile  G(r/) ,  the  resulting  velocity  field  is  found  as: 
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u(ri())=m=m m 

r  r 


Also  the  pressure  field  can  be  found  as  follows: 

dp  _2//F(0)G'  dpjrpj)  _  2//F(0)G' 
SO  r 2  dij  r 2 


(20) 


(21) 


Integrating  once  with  respect  to  6  or/; , 


/?(r,  /;)  =  G(/;)  +  c(r) 


Differentiating  with  respect  to  r, 


|L=zf^)G(,)+cXr) 

or  r 


(22) 


The  second  momentum  equation  solved  for  the  pressure  gradient  gives: 


dp  _pF2  |  //Fff  ^  d//(r,/;)_pF(0)2G2  |  //F(0)Gff 


dr  r3  r3 


dr 


2  3 

a  r 


(23) 


Equating  the  above  last  two  equations  and  solving  for  the  constant  c(r)  we  get: 


-4//F(0) 


G(/7)  +  c'(r)  = 


pF(0)2G 2  pF(0)G" 


+  - 


a2r3 


(24) 


c\r)  = 


a  r 
pF(  0) 


2  3 

a  r 


G"+a  F(Q)G2+4a2G 


G"  +  a  Re  G2  +  4a2G 

G^O) 


//F(0)  , 

c(r)=-v^G  (1)  +  c> 

2a  r 


(25) 


Thus  the  pressure  field  is  given  by: 
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r  2a  r 


p(r,V>  Wg^-G'XI)  ]  +  c, 

/  rt  !/•  L  -l 


2a  r 


(26) 


But, 


pFj 0)  _  /£  2  _  V7  2  2  1 

2a2r2  _  2ZF(0)  Mm  _  2a2F(0)  _  2  aRe 


And  if 


r  — » oo  =>  p(r,T])^>cl=px, 


This  means  pressure  at  infinity  becomes  uniform  i.e.  the  diffuser  is  very  long. 

Thus, 


P(r,v)~ px=^pu2m  — [4a2G(/7)-G"(l)  ]  => 


Poo  -p(r,tj)  i 


a  Re 


[G'(l)-4a2G(7)  ] 


(27) 


This  is  the  pressure  coefficient  that  specifies  the  recovery  in  static  pressure  as  a 
function  ofr,p. 

2.  Solution  of  the  Problem 


The  nonlinear  equation  (15)  will  be  solved  in  MATLAB  using  the  build  in 
RUNGE-KUTTA  45  algorithm  for  ordinary  differential  equations.  Because  of  the 
accuracy  and  robustness  of  this  scheme  the  solution  can  be  regarded  close  to  an  exact 
analytical  solution.  The  programming  specifics  are  given  in  Appendix.A.  Here  the  results 
are  presented  in  the  form  of  velocity  profiles  for  both  nozzles  and  diffusers. 
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Figure  8.  Velocity  profiles  for  various  Reynolds  numbers  in  converging,  parallel  and 

diverging  channels. 

It  is  seen  in  Figure  8  that  multiple  solutions  are  possible  depending  on  the 
Reynolds  number.  Both  symmetric  and  asymmetric  velocity  profiles  are  possible  in  a 
diffuser.  The  elliptic  integral  solution  [21]  of  the  non-linear  differential  equation  (15) 
predicts  a  separation  criterion  of  a  Re  >10.31.  For  a  five  degree  half  angle  divergent 
walls  this  means  that  separation  is  possible  if  Re  >  118. 

C.  NUMERICAL  FORMULATION  AND  SOLUTION  OF  THE  PROBLEM  IN 
CFD-ACE+. 

1.  Formulation  of  the  Problem 

The  software  CFD-ACE  must  be  supplied  with  a  specific  bounded  geometry. 
Here  we  choose  a  diffuser  of  half  angle  5  degrees  with  the  spatial  dimensions  shown  in 
the  following  figure  in  [m]: 
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The  area  ratio  is  approximately  3  and  the  non-dimensional  length  23.  This  gives  a 
five-degree  half  angle.  The  objective  is  to  compute  the  flow  field  in  CFD-ACE  for 
boundary  conditions  taken  the  same  as  the  ones  used  in  the  analytical  solution.  The 
computed  flow  field  will  be  checked  against  the  analytical  field  at  two  stations  for  the 
velocity  profiles  (0.004  m  and  0.007  m)  and  the  centerline  for  the  static  pressure 
distribution.  A  300X120  grid  and  air  as  a  medium  will  be  used  in  the  numerical  solution. 
The  Reynolds  number  will  be  chosen  to  be  Re=  100/0.08726=1 146,  a  value  which  gives  a 
dimensionless  profile  approaching  separation  as  shown  in  Figure  8. 

The  dimensional  velocity  and  pressure  profiles  needed  for  boundary  conditions 
are  formulated  in  MATLAB  and  for  the  inlet  are  given  graphically  in  Figure  10.  From  the 
definition  of  the  Re  and  the  velocity  profiles  at  the  inlet  and  outlet  of  the  diffuser  we  get: 


Figure  10.  Velocity  profiles  at  the  inlet.  On  the  left  non-dimensional,  right  top 
dimensional  and  right  middle  and  bottom  decomposed  to  dimensional  x  and  y 

components  [m/sec] . 
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2.  Solution  of  the  Problem 

The  solution  of  ACE+  is  visualized  in  CFD-VIEW.  Results  are  presented  in  the 
forms  of  velocity  contours  streamline  contours  and  static  pressure  contours,  in  Fig.9,  10, 
1 1  respectively. 


Figure  13.  Static  pressure  contours. 

Comparison  of  the  numerical  to  the  analytical  solution  in  terms  of  velocity 
profiles  is  given  in  Fig.9.  This  plot  presents  the  non-dimensional  velocity  profiles  at  two 
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stations  namely  at  x=0.004  m  and  x=0.007  m  along  with  the  exact  velocity  profile  for  the 
same  mass  flow  rate.  The  agreement  between  the  two  solutions  is  excellent  with 

Au/  <0.1%. 

/U 


0/  a 


Figure  14.  Comparison  of  the  computed  velocity  profiles  with  the  exact  solution 


Comparison  of  the  numerical  to  the  analytical  solution  in  terms  of  static  pressure 
profiles  along  the  centerline  is  given  in  Fig.  15.  The  difference  8P is  from  a  value  of 
101,325  Pa  assumed  at  the  diffuser  exit.  Again  the  agreement  between  the  two  solutions 
is  satisfactory  (~  1%). 
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Figure  15.  Comparison  of  the  computed  static  pressure  centerline  profiles  with  the  exact 

analytical  solution 


D.  A  NOTE  ON  THE  NON-EXISTENCE  OF  SIMILARITY  SOLUTIONS  IN 
CONICAL  DIFFUSERS  WITH  POINT  SOURCE  ENTRY  FLOWS 
SINGLE  SPACE 

The  possibility  of  existence  of  similarity  solutions  in  conical  diffusers  is 
investigated  here.  If  the  flow  is  assumed  purely  radial  emerging  from  a  point  source  at  the 
origin  and  driven  through  a  conical  diffuser  the  analogous  3-D  case  to  the  Jeffery-Hamel 
flows  in  two  dimensions  emerges.  Equations  of  continuity  and  momentum  are  written  in 
the  spherical  coordinate  system: 
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(28) 


1  d(r2ur)  |  1  (ug  sing)  |  1  du+  _Q 

r2  dr  rsinQ  d6  rsin#  d(j) 


DV  „  — 

P ~Dt=~Vp  +  ^  V  (29) 


Since  the  velocity  emerging  from  the  point  source  is  assumed  to  be  purely  radial, 
the  ug  ,  uifi  components  of  the  velocity  must  be  zero.  Substituting  in  continuity  above  and 

integrating  once  with  respect  to  r  we  get: 


WW) 


(30) 


Ug=0  (31) 

«+=  0  (32) 


We  search  for  a  similarity  velocity  profile  F(0,  (/>) .  It  will  be  shown  in  two  ways 
that  such  a  solution  does  not  exist.  First  by  dimensional  reasoning  (presented  here)  and 
then  in  a  more  mathematical  form  presented  in  the  Appendix  B.  The  velocity  at  the 
centerline  is 


um=ur(rM) 


vF(  0,0 


A  Reynolds  number  can  be  defined  as, 

Rc  _  L  rl l r  -  F(°>  a 

v  v  r 


Where  a  is  the  diffuser  half  angle. 

This  number  gives  the  ratio  of  inertia  to  viscous  forces  inside  the  diffuser.  If 
F(0,  </))  were  a  similarity  profile,  then  F(0,  </>)  at  the  centerline  would  be  invariant.  Since 
a  is  constant,  Re  number  is  a  function  of  r.  Thus  Re  number  is  not  invariant  meaning,  the 
ratio  of  the  inertia  to  viscous  forces  decreases  downstream.  Thus  similarity  solutions  to 
the  Navier-Stokes  equations  for  this  type  of  flows  do  not  exist.  On  the  contrary  this  is  not 
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the  case  for  the  2-D  plane  diffusers  presented  in  the  previously.  There,  the  Reynolds 
number  is  defined  as, 


Re  = 


ura 


v 


(33) 


And  the  centerline  velocity  as, 


u 


m 


F{0) 


(34) 


Substituting  in  (21)  above, 


F{0) 

— T7 1  n\ 

r  _  F(9)a 
v  v 


(35) 


Thus  Reynolds  number  is  clearly  invariant  throughout  the  flow  field  since  a  and 
v  are  constants.  This  means  similarity  solutions  are  possible  in  this  case. 
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IV.  NUMERICAL  SIMULATION  OF  INCOMPRESSIBLE 
LAMINAR  STEADY  FLOW  IN  2-D  DIFFUSERS  WITH  EXIT 

CHANNELS 


A.  OVERVIEW 

Direct  numerical  simulation  will  be  presented  to  investigate  performance  and 
various  flow  regimes  that  might  exist  in  2-D  plane  laminar  diffuser  flows  at  Reynolds 
numbers  range  of  105-1048.  Channels  are  situated  upstream  and  downstream  of  the 
diffuser.  A  number  of  different  geometries  (AR=1.15  through  5  and  L/Wi=l-48)  are 
considered.  The  diffuser  inlet  flow  is  assumed  incompressible,  laminar  and  fully 
developed. 

B.  DEFINITION  OF  THE  PROBLEM 

1.  Geometric  Configuration 

The  general  configuration  used  in  this  study  and  the  appropriate  geometric 
parameters  that  describe  it  are  shown  in  Fig.  17.  In  this  figure  Wi,  W2  are  the  inlet  and 
outlet  heights  of  the  diffuser,  L  is  the  length  of  the  diffuser. 


The  area  ratio  is  varied  from  1 . 15  to  5  and  the  length  to  inlet  height  L/Wi.is  varied 
from  1  to  48.  Inlet  and  outlet  channel  lengths  are  chosen  with  appropriately  large  length 
to  accommodate  properly  the  flow  inlet  and  outlet  boundary  conditions.  It  is  explained 
later  why  and  how  these  lengths  are  chosen. 
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2.  Simulations  Matrix 

For  the  geometries  mentioned  in  1  above,  Reynolds  number  runs  of  63,  84,  105, 

210,  314,  420,  629,  839  and  1048  will  be  performed.  The  Reynolds  number  is 
defined  as: 

Re  =  ^ 

V 


Where  U  is  the  mean  inlet  channel  velocity,  Wi  is  the  inlet  channel  height  and  v 
is  the  kinematic  viscosity.  The  specific  test  geometries  in  terms  of  area  ratio  and  non- 
dimensional  length  are  shown  in  Figure  18  on  a  log- log  graph. 


The  same  information  in  different  form  is  given  in  Table.l  where  the  diffuser  half 
angles  in  degrees  are  included.  Geometries  of  half  angles  more  than  45°  and  less  than  1° 
are  discarded.  At  the  45°  angles  the  grids  become  skew  enough  to  influence  the  accuracy 
and  convergence  of  the  solution  and  at  1°  and  below  geometries  are  of  no  practical 
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interest.  For  all  diffuser  geometries  the  flow  is  computed  for  the  different  Re  numbers 
mentioned  before.  Every  such  run  will  be  referred  to  as  a  “case”  from  now  on.  For 
instance  a  diffuser  of  non-dimensional  length  4.5,  area  ratio  3.5  and  Re  number  314  is  the 
“case”  4.5  3.56.  Table. 2  contains  Re  numbers  and  max  centerline  velocities  for  all  cases. 


Table  1  Diffuser  geometries  tested. 


AR— » 

1.1 

1.15 

1.25 

1.5 

2 

2.5 

3 

3.5 

4 

4.5 

5 

Set48 

1.19 

1.49 

1.79 

2.08 

2.38 

Set36 

1.19 

1.59 

1.98 

2.38 

2.78 

3.18 

Set24 

1.19 

1.79 

2.38 

2.98 

3.57 

4.17 

4.76 

Setl8 

1.59 

2.38 

3.18 

3.97 

4.76 

5.55 

6.34 

Setl2 

1.19 

2.38 

3.57 

4.76 

5.95 

7.12 

8.3 

9.46 

Set9 

1.59 

3.18 

4.76 

6.34 

7.91 

9.46 

11 

12.53 

Set6 

1.19 

2.38 

4.76 

7.12 

9.46 

11.7 

14.03 

16.3 

18.4 

Set4.5 

1.59 

3.18 

6.34 

9.46 

12.5 

15.5 

18.4 

21.5 

23.9 

Set3 

1.43 

2.38 

4.76 

9.46 

14.03 

18.4 

22.62 

26.5 

30.2 

33.7 

Setl.5 

1.91 

2.86 

4.76 

9.46 

18.43 

26.56 

33.69 

39.8 

45 

Setl 

2.86 

4.29 

7.12 

14.03 

26.5 

36.87 

45 

Table  2  Reynolds  numbers  tested  for  each  diffuser  geometry. 


Case# 

1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

Re 

21 

63 

84 

105 

210 

314 

420 

629 

839 

1048 

Umax 

(m/sec) 

0.5 

1.5 

2 

2.5 

5 

7.5 

10 

15 

20 

25 

3.  Grid  Types 

Grids  are  generated  for  each  case  in  CFD-GEOM  and  then  imported  in  the  solver. 
The  two  general  arrangements  are  given  schematically  in  Fig.  19,  20: 
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F igure  1 9 .  Grid  type  1 . 


Figure  20.  Grid  type2. 


Type  2  grids  are  used  whenever  the  diffuser  exit  comer  angles  of  type  1  grids 
become  less  than  70°  i.e.  for  geometries  where  half  diffuser  angle  exceeds  20°. 

4.  Boundary  Conditions 

The  types  of  boundary  conditions  applied  here  are: 

(1)  The  no  slip,  no  penetration  condition  at  the  walls. 

(2)  A  fully  developed  parabolic  velocity  profile  at  the  inlet  of  the  upstream 
channel  calculated  in  MATLAB  and  imported  in  ACE+.  The  max  centerline  velocity  of 
this  profile  is  varied  according  to  Table.2. 
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(3)  A  uniform  static  pressure  condition  at  the  outlet  of  the  downstream  channel. 

To  determine  a  good  location  for  the  boundaries  of  the  computational  domain 
inlet  and  outlet  so  that  the  boundary  conditions  would  be  satisfied  sufficiently,  the 
following  guidelines  were  followed.  For  the  upstream  boundary,  where  the  parabolic 
profile  would  be  imposed,  the  streamlines  should  have  to  be  parallel.  Separate 
computations  (performed  in  MATLAB  using  a  finite  difference  method)  showed  that  for 
potential  flow  through  a  sudden  expansion  with  area  ratio  5  (the  max  area  ratio  used  in 
this  study)  the  streamlines  were  parallel  at  a  distance  5Wi  upstream  the  expansion.  Thus 
an  upstream  channel  length  of  8W)  was  chosen  throughout  this  work.  For  the 
downstream  exit  boundary  the  streamlines  should  again  be  parallel  for  the  uniform 
pressure  to  be  imposed.  Thus  channel  lengths  of  80Wi-120Wi  are  used,  based  on 
preliminary  studies,  done  by  the  writer  and  several  other  studies  such  as  in  ref.  [18].  See 
Appendix.C  for  the  MATLAB  coding  and  results  on  the  potential  flow  through  the 
sudden  expansion. 

C.  NUMERICAL  SOLUTION  IN  ACE+ 

1.  Description  of  Solver  Settings 

The  “flow”  module  of  ACE+  is  used.  That  means  continuity  and  N-S  equations 
will  be  solved  numerically  in  two  dimensional,  steady,  incompressible  forms.  Air  at 
atmospheric  conditions  will  be  used  as  the  medium.  Inlet  channel  max  velocities  will  be 
of  the  order  of  1.5  to  20  m/sec.  To  keep  flow  laminar  and  Re  in  the  desired  interval,  inlet 
channel  width  is  chosen  to  be  1  mm.  All  other  lengths  are  scaled  accordingly  in 
millimeters.  Guess  values  for  velocity  and  pressure  should  be  used  to  initialize  the  flow 
field.  Following  the  manual,  the  exit  static  pressure  is  used  as  initial  value  for  the 
pressure  field  in  order  to  achieve  faster  convergence.  A  value  for  velocity  from  Table.2  is 
used  throughout  depending  on  the  particular  run. 

Numerical  results  are  digitally  stored  in  DTF  files,  named  as  follows:  A 
diffuser  of  non-dimensional  length  4.5,  area  ratio  3.5and  Re  number  320  is  the  file 
4.5  3.56.DTF.  The  last  digit  specifies  Re  number  from  Table.2.  The  folder  containing  all 
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the  runs  for  the  4.5  3. 5  geometry  specified  above  would  be  named  D4.53.5.  Thus  this 
folder  would  contain  DTF  files  such  as  4.5  3.54.DTF,  4.5  3.55.DTF, 

4.5_3.56.DTF. .  .etc. 

2.  Convergence  of  the  Numerical  Solution 

As  already  mentioned  above,  a  number  of  runs  will  be  performed  for  the  various 
geometries  and  Reynolds  numbers.  For  each  “case”  the  residual  history  of  the  dependent 
variables  and  the  “mass  flow  summary”  was  recorded.  Depending  on  diffuser  geometry 
and  Re  number,  a  general  trend  was  singled  out  throughout  the  computations. 
Specifically,  for  certain  combinations  of  diffuser  geometry  and  Re  number  the  residual 
curves  obtained  had  a  form  shown  in  the  following  Figure  21 : 


In  these  cases  the  computer  programm  neither  converged  nor  diverged.  The  “mass 
flow  summary”  was  very  poor  meaning  continuity  was  satisfied  no  better  than  two  to 
three  orders  of  magnitude  and  generally  more  that  six  orders  is  sought.  Velocity  contours 
for  two  successive  iterations  of  the  same  case  are  shown  in  Figure.22.  It  is  evident  that 
the  flow  field  changes  back  and  forth  without  converging  to  a  particular  form. 
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Figure  22.  Non-converging  case  velocity  contours. 


For  other  cases  the  solution  converges  in  less  than  1500  iterations  and  the 
residuals  have  droped  at  least  seven  orders  of  magnitude.  Mass  flow  summaries  are  very 
good  indicating  continuity  satisfaction  of  at  least  nine  to  ten  orders  of  magnitude.  Figures 
23  and  24  show  residuals  and  u  velocity  contours  for  such  cases  respectively: 


Figure  23.  Typical  residual  history  in  steady  converging  cases. 
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Figure  24.  Typical  converging  case  velocity  contours.  Top  symmetric  separated  flow. 

Bottom  asymmetric  separated  flow. 


The  distinct  trend  of  the  residuals  and  mass  flow  summaries  is  observed  in  all 
steady  computations.  As  will  be  discussed  later  this  trend  will  be  the  decisive  factor  for 
distinguishing  flow  regimes. 

3.  Grid  Independence  Study 

All  numerical  calculations  are  performed  solving  the  descretized  form  of  the 
appropriate  differential  equations  on  a  finite  grid.  It  is  well  known  that  the  density  of  the 
grid  has  a  significant  impact  on  the  accuracy  and  convergence  of  the  computations.  Thus 
we  need  to  perform  a  grid  independence  study  until  a  grid  density  is  found  that  gives  a 
grid  independent  solution. 

In  this  study  computations  on  various  grid  densities  were  carried  out  at  Re 
numbers  420  and  630  on  diffusers  of  non-dimensional  length  L/Wi=  24  and  area  ratios  of 
3.5,  3  and  2.5,  3  respectively.  Two  diffuser  grids  are  examined.  One  more  dense  and  one 
less  dense  than  the  nominal  one  used  for  the  bulk  of  the  computations.  Since  the  nominal 
24set  diffuser  grid  dimensions  are  (500X80),  a  denser  grid  (750x100)  and  a  coarser  grid 
(300x60)  are  employed  for  grid  independence  study.  The  results  are  compared  to  see  if 
the  solutions  deviate  from  one  another.  It  is  found  that  good  agreement  between  the  three 
solutions  exists.  Velocity  profiles  and  centerline  static  pressure  distributions  will  be 
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compared  for  cases  where  convergence  is  achieved.  For  cases  where  convergence  is  not 
achieved  the  residual  histories  will  be  compared. 

The  results  of  the  grid  independence  study  are  shown  in  Figures  25  and  26.  Good 
agreement  is  observed  between  the  solutions  (relative  difference  less  than  1%) 


u/u 


max 


Figure  25. 


Velocity  profiles  comparisons  with  different  grid  density. 
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Figure  26.  Comparison  of  centerline  static  pressure  coefficient  in  terms  of  non- 
dimensional  length  for  grids  of  different  spatial  resolution. 


Figure  27.  Residuals  vs.  iteration  number  comparison  for  the  coarse  and  dense  grid  non¬ 
converging  cases  at  Re=629  respectively. 
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D.  COMPUTATIONAL  RESULTS  AND  DISCUSION 


Residual  histories  of  the  primary  dependent  variables  (two  velocity  components 
and  pressure)  and  the  “mass  flow  summary”  are  the  criteria  distinguishing  flow  regimes. 
Cases  where  residuals  are  of  the  form  of  Figures  21  and  27  and  mass  flow  summary  is 
poor  are  characterized  “non  2-D”,  meaning  that  the  flow  is  no  longer  two-dimensional  or 
time  independent.  Although  the  boundary  conditions  are  time  independent,  the  flow  is 
unsteady  and/or  3-D  and  as  a  result  the  steady  two-dimensional  equations  are  no  longer 
suitable  to  describe  the  phenomena  inside  and  downstream  the  diffuser.  Cases  where  the 
residuals  are  of  the  form  of  Figure  23  are  characterized  steady.  The  flow  is  two- 
dimensional  steady,  separated  from  one  or  both  diffuser  walls  except  for  very  low  angle 
geometries  were  it  is  not  separated  at  all. 

1.  Flow  Regime  Maps 

Figure  28  depicts  a  simple  sketch  of  the  different  flow  regimes  identified  in  this 
study  for  a  single  Re  number. 


Flow  regime  map  sketch 


Not  converging 


Figure  28.  Flow  regime  map  sketch  for  Reynolds  number  of  3 14. 
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As  the  divergence  angle  is  increased  there  are  four  regions  of  different  flow 
configurations  observed.  Starting  at  the  arrow  for  example  at  L/Wi=12  the  first  solid  line 
marks  the  attached  steady  flow  regime  A.  Increasing  the  angle  of  divergence  we  enter 
region  B  were  the  flow  in  the  diffuser  is  symmetrically  separated.  Increasing  the  diffusion 
angle  further,  region  C  is  reached,  where  the  flow  is  separated  but  asymmetric  relative  to 
the  diffuser  centerline.  Above  the  line  between  C  and  D  is  the  region  where  the  numerical 
solution  of  the  2-D  equations  neither  converges  nor  diverges.  The  lines  between  A  B  and 
C  regions  are  qualitatively  drawn  that  far  apart  because  attached  and  symmetric  separated 
flows  occur  for  very  low  divergence  angles  (2-4  degrees)  so  in  a  quantitative  graph  they 
would  be  very  close.  Region  C  is  the  one  that  takes  up  the  largest  space  in  the  map  below 
the  Re=314  line. 

The  line  between  regimes  C  and  D  is  the  one  computed  quantitatively  and  is  given 
in  Figure  29  for  six  different  Reynolds  numbers.  Appendix  D  explains  the  process  by 
which  Fig.29  is  extracted  from  the  numerical  data.  Above  the  lines  in  Fig.29,  where  there 
is  no  convergent  2-D  steady  solution,  the  flow  is  believed  to  be  three-dimensional  and 
unsteady.  Figure  22  depicts  velocity  contours  of  such  flow  where  it  can  be  seen  that  the 
contours  form  an  unphysical  configuration  that  more  likely  depicts  a  snapshot  of  three- 
dimensional  unsteady  flow  vertical  center  plane  cut. 
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Figure  29.  AR-L/Wi  flow  regimes  map. 


A  different  form  of  Fig.29  is  given  in  Fig.30.  In  this  figure  the  power  laws  derived 
for  the  lines  in  Fig.29  and  presented  in  Appendix  D  are  used  to  correlate  the  data  in  a  plot 
of  AR  vs.  Re  number  with  parameter  L/Wi.  Notice  that  the  lines  out  of  the  box  in  Fig.30 
are  extrapolations  to  the  data  and  should  be  used  with  caution. 

Contrary  to  the  high  Re  number  regime  map  which  only  shows  geometry 
dependence  of  the  flow  regimes,  at  low  Re  numbers  with  fully  developed  laminar  flows 
entering  the  diffuser,  the  flow  regimes  depend  strongly  on  diffuser  geometry  and 
Reynolds  number.  With  increasing  angle  the  flow  goes  from  attached  to  symmetric 
separated  to  asymmetric  separated  and  finally  to  a  configuration  that  is  not  2-D  but  is 
believed  to  be  three  dimensional  and/or  unsteady.  There  is  no  region  of  transitory  stall 

present  between  the  attached  and  2-D  asymmetric  separated  regimes. 
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Figure  30.  AR-Re  flow  regimes  map. 


2.  Flow  Visualization  Graphs 

To  gain  a  better  visual  understanding  on  the  different  flow  configurations 
described  on  the  regime  map  the  following  graphs  may  be  used.  Figures  31,  32,  33  and 
34  depict  velocity,  streamline  and  static  pressure  contours  respectively.  Point  tracing 
emphasizes  recirculation  regions  of  the  flow.  For  low  diffusion  angles  the  flow  is 
attached,  then  goes  to  separated  symmetric,  then  to  separated  asymmetric  and  finally  it 
takes  the  configuration  of  Fig.33  where  the  2-D  flow  brakes  down  to  what  is  believed  to 
be  the  onset  of  three  dimensional  and/or  unsteady  flow.  The  pressure  field  is  for  the 
symmetric  cases  is  fairly  uniform  but  becomes  less  uniform  for  the  asymmetric  separated 
case. 
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Figure  3 1 .  Velocity,  streamlines  and  static  pressure  contours  for  attached  flow. 


Figure  32.  Velocity,  streamlines  and  static  pressure  contours  for  symmetric  separated 

flow. 
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Figure  33. 


Figure  34. 


Velocity,  streamline  and  static  pressure  contours  for  asymmetric  separated 

flow. 


Velocity,  streamline  and  static  pressure  contours  for  a  non-converging  case. 
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3.  Pressure  Recovery  Maps 

A  pressure  coefficient  Cp  defined  in  equation  (22)  is  used  to  evaluate  pressure 
recovery  along  the  diffuser-exit  channel  centerline  in  the  flow  direction: 


Cp(x)  = 


/?(*)- A 

1/2  pU\ 


(36) 


The  static  pressures  in  all  steady  2-D  cases  (both  symmetric  and  asymmetric 
flows)  were  found  to  be  uniform  or  almost  uniform  in  the  y  direction.  Thus  centerline 
pressure  distribution  is  taken  to  calculate  the  value  of  pressure  recovery.  A  typical  Cp  vs. 

xtWx  graph  is  shown  in  figure.31.  In  most  of  the  cases  tested  the  flow  is  separated  from 
the  diffuser  walls  and  reattaches  in  the  downstream  channel.  These  reattachment  lengths 
are  twice  three  times  or  even  longer  than  the  diffuser  length.  Thus  the  Cp  chosen 

represents  performance  of  the  diffuser  and  downstream  channel  together. 


Figure  35.  Typical  pressure  recovery  for  different  Re  numbers  along  the  diffuser-exit 

channel  centerline. 
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Several  features  in  these  charts  should  be  pointed  out  (also  discussed  in  [18]). 
First  the  point  of  maximum  Cp  does  not  occur  at  the  diffuser  exit  but  in  the  down  stream 

channel  near  the  reattachment  point.  Maximum  Cp  always  shifts  to  the  left  i.e.  further 

downstream  with  increasing  Re  number  since  the  reattachment  point  does  so  also.  The 
Cp  maximum  does  not  always  increase  as  the  Re  number  increases  while  keeping  the 

same  geometry  but  reaches  a  max  value  at  intermediate  Re  numbers.  This  implies  the 
existence  of  an  optimum  performance.  The  wavy  pattern  of  the  curves  and  the  consequent 
overall  drop  of  recovery  is  the  result  of  the  asymmetrical  flows.  As  noted  earlier  such 
flow  conditions  are  observed  approaching  the  regime  curves  from  below. 

Diffuser-exit  channel  operating  optimums  can  be  extracted  from  the  pressure 
recovery  curves.  Such  performance  maps  are  constructed  if  the  maximum  recovery  Cp  m 

is  plotted  against  area  ratio  and  non-dimensional  length  with  Re  numbers  as  parameter 
and  after  that  the  AR-L/W1  cross-plots  are  taken.  Maps  at  Re=105,  210,  315,  420  and 
629  were  created.  For  Re=105,  629  the  maps  are  presented  here.  The  rest  three  with  their 
respective  Cp-AR,  Cp-L/Wi  graphs  can  be  found  in  Appendix  E.  For  Re  numbers  above 
629  not  enough  numerical  data  exist  because  of  the  limitations  imposed  by  the  respective 
flow  regime  curves. 

In  Figure  38  and  41  the  dotted  lines  represent  the  respective  flow  regime  curves 
on  the  performance  maps.  The  optimum  diffuser  performance  exists  always  below  these 
lines.  Thus  in  the  case  of  low  Re  numbers  the  max  diffuser  performance  does  not  exist  in 
the  unsteady  regime  which  is  above  the  curves  but  below  the  regime  lines  in  the  two 
dimensional  asymmetric  flow  regime.  The  opposite  happens  in  high  Re  flows  where  a 
diffuser  in  the  transitory  stall  regime  exhibits  the  best  pressure  recovery. 

The  flow  regime  lines  in  Figure  29  are  approximately  lines  of  constant  diffusion 
angle.  In  this  context  maximum  pressure  recovery  in  low  Reynolds  number  diffusers 
takes  place  for  increasing  diffusion  angle  as  Re  decreases.  This  is  in  contrast  to  high  Re 
diffusers  where  the  maximum  recovery  comes  for  a  single  diffusion  half  angle  of  6-7 
degrees. 
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Figure  36.  Variation  of  maximum  Cp  with  AR  for  Re=105. 


Figure  37.  Variation  of  maximum  Cp  with  L/Wi  for  Re=105 
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F igure  3  8 .  AR-L/W i  performance  map  for  Re=  105. 


Figure  39.  Variation  of  maximum  Cp  with  L/Wi  for  Re=629. 
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Figure  40.  Variation  of  maximum  Cp  with  AR  for  Re=630 


Figure  41.  AR-L/Wi  performance  map  for  Re=629. 
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As  Reynolds  number  increases  the  optimum  region  shifts  to  lower  diffusion 
angles  and  higher  diffuser  lengths.  As  diffuser  lengths  become  smaller  recovery  should 
go  to  a  value  close  to  that  of  an  abrupt  expansion  with  the  same  area  ratio.  In  Figures  37 
and  39  the  pressure  recovery  on  the  left  approaches  a  single  value  for  each  area  ratio.  A 
simple  pseudo-inviscid  analysis  presented  in  Appendix  C  suggests  that  the  maximum 
pressure  recovery  achieved  with  an  abrupt  expansion  is  0.5  when  the  area  ratio  is  2.  Thus 
Figures  37  and  39  suggest  that  it  is  advantageous  to  use  diffusing  passages  instead  of 
expansions  in  order  to  get  higher  pressure  recovery  when  the  Reynolds  number  is  low. 
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V.  TIME  RESOLVED  SIMULATION  OF  2-D  FLOWS  IN 
DIFFUSERS  WITH  EXIT  CHANNELS 

A.  OVERVIEW  AND  PURPOSE 

A  time  marching  time  accurate  numerical  solution  of  the  unsteady  2-D  Navier- 
Stokes  equations  was  attempted  here.  The  flow  was  impulsively  started  with  a  pressure 
force  and  after  sufficient  time  approaches  a  condition  that  depends  on  diffuser  geometry 
and  Re  number.  It  is  attempted  to  see  if  the  flow  approaches  the  same  2-D  condition  as 
when  solving  the  steady  equations,  for  AR-L/Wi-RE  combinations  situated  upon  and 
below  the  two-dimensional  flow  regime  curves  computed  in  the  previous  chapter.  The 
same  time  marching  procedure  is  applied  to  AR-L/Wi-RE  combinations  that  lay  above 
the  regime  lines  in  an  attempt  to  see  what  happens  to  the  flow  field  after  a  large  period  of 
time. 

B.  DEFINITION  OF  THE  PROBLEM 

1.  Geometry  Configuration,  Simulations  Matrix  and  Grid  Types 

The  same  basic  Fig.  15  geometry  is  considered  again  for  the  time  resolved 
solutions.  Inlet  channel  length  is  modified  for  reasons  explained  later.  The  grid  types  are 
the  same  as  the  ones  used  in  the  steady  computations.  Different  runs  will  be  performed 
according  for  the  cases  described  in  Table.3.  A  letter  (t)  denotes  a  transient  case.  The  rest 
is  being  kept  the  same  as  the  steady  case  notation. 


Table  3  Time  resolved  cases  simulations. 


Re=314 

3_36t 

4.5_3.56t 

6_3.56t 

9_46t 

Above  b-b 

Re=314 

3_2.56t 

4.5_36t 

6_36t 

9_36t 

Below  b-b 

Re=210 

4.5_55t 

1.5_45t 

Above  a-a 

Re=210 

4.5_4.55t 

1.5_3.55t 

Below  a-a 

On  the  flow  regime  maps  the  test  geometries  show  as  in  Fig.41.  In  this  figure  the 
triangles  depict  cases  given  in  Table.3  above  and  below  the  flow  regime  lines  for 
Reynolds  numbers  210,  314. 
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Figure  42.  Unsteady  cases  on  the  flow  regime  map. 

2.  Boundary  and  Initial  Conditions 

For  the  time  resolved  simulations  a  constant  total  pressure  is  used  at  the  inlet 
boundary  of  the  upstream  channel.  At  the  downstream  channel  exit  a  constant  static 
pressure  boundary  is  imposed.  An  initial  condition  of  static  pressure  equal  to  the  exit 
static  pressure  and  zero  velocity  throughout  is  used.  The  inlet  channel  length  needs  to  be 
modified  to  allow  fully  developed  parabolic  flow  at  the  diffuser  entrance.  An 
experimentally  determined  correlation  of  the  following  type  (see  Ref.3)  is  applied  to 
determine  the  appropriate  upstream  channel  length: 

x  0.63 

— = - +  0.044  Re 

h  0.035  Re+1 

Thus  for  Reynolds  number  210,  314  minimum  lengths  of  x  /  /? « 1 0  and 
x/A  « 14  are  used  respectively. 
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c. 


NUMERICAL  SOLUTION  IN  ACE+ 


1.  Description  of  ACE  Settings  for  the  Transient  Computations 

The  flow  module  is  again  used  for  the  transient  simulations.  Specifying  the  time 
step  and  total  number  of  steps  in  the  appropriate  fields  the  unsteady  two-dimensional 
incompressible  Navier-Stokes  equations  are  solved.  For  all  calculations  a  time  step  of 
0.0015  sec  is  used.  Maximum  number  of  time  steps  is  computed  by  the  max  time  the  flow 
needs  to  approach  the  steady  condition.  An  estimate  of  max  time  needed  by  the  flow  to 
approach  the  steady  flow  can  be  taken  by  an  existing  exact  Navier-Stokes  solution  for  the 
impulsively  starting  flows  in  channels  due  to  pressure  gradients.  This  exact  solution  gives 
the  following  non-dimensional  time  parameter: 


T 


tv 


h2 


For  channel  flows  it  is  found  that  the  velocity  approaches  the  parabolic  shape  at: 

r  -0.47 


Thus  for  our  case,  an  estimate  for  the  max  time  the  flow  needs  to  approach  a  final 
steady  condition  might  be: 


0.47  h2 

t= - 

v 


Calculations  herein  use  the  viscosity  of  air  1.589e-5  m2/sec  and  the  height  of  the 
exit  channel.  A  15-20%  increase  to  the  above  calculated  time  determines  the  actual  value 
that  is  finally  used.  The  overall  geometry  dimensions  are  in  millimeters.  This  ensures  that 
the  pressure  gradient  is  felt  instantly  throughout  the  flow  field  because  the  time  step  used 
is  at  least  three  times  more  than  the  time  sound  waves  need  to  propagate  downstream. 
Under  relaxation  values  of  0.06-0.08  are  used  for  all  dependent  variables. 

2.  Convergence  of  the  Transient  Solutions 

Convergence  history  of  the  dependent  variables  in  each  time  step  can  be  viewed 

in  ACE+.  For  cases  below  the  lines  a-a  or  b-b  on  Figure  42,  residuals  steady  out  after 

sufficient  amount  of  time.  For  cases  above  the  lines  this  is  not  true.  Residuals  continue  to 

converge  to  the  specified  criterion  in  each  time  step  and  the  flow  field  solution  seems  to 

alternate  between  several  configurations  (in  each  time  step).  This  is  believed  to  be  the 

51 


onset  of  three-dimensional  and/or  unsteady  flow.  The  following  residual  curves  are 
obtained  for  large  periods  of  time  depending  on  whether  the  geometry  is  above  or  below 
the  flow  regime  curve  of  Figure  42. 


Figure  43.  Typical  time  step  residual  history  at  large  times  for  cases  below  a-a  or  b-b. 


Figure  44.  Typical  time  step  residual  history  at  large  times  for  cases  above  a-a  or  b-b. 
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D. 


COMPUTATIONAL  RESULTS  AND  FLOW  VISUALIZATION. 


Transient  computations  results  will  be  presented  in  terms  of  evolution  in  time  of 
velocity  profiles  and  centerline  static  pressure  distributions.  These  results  will  be 
compared  with  analogous  ones  computed  in  the  previous  chapter  for  the  same  geometries 
and  Re  numbers  only  for  the  steady  flow  cases.  Velocity  and  streamline  contours  will  be 
presented  for  all  cases  (both  above  and  below  the  lines  (a-a)  and  (b-b). 

1.  Instantaneous  Velocity  and  Centerline  Pressure  Profiles 
Velocity  distributions  in  the  upstream  and  far  downstream  in  he  exit  channel 
approach  the  parabolic  profiles  when  the  steady  flow  is  reached.  The  following  Figure  44 
depicts  a  typical  evolution  of  velocity  profile  at  the  exit  of  the  downstream  channel. 


u/u 


Figure  45.  Instantaneous  and  final  steady  velocity  profiles. 

The  non-dimensional  time  constant  for  the  flow  to  approach  the  parabolic  profile 

is  found  to  be  0.1 1  for  downstream  channel  height  corresponding  to  a  diffuser  of  AR=3. 
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The  circles  depict  the  exit  velocity  distribution  for  the  equivalent  case  of  same  Reynolds 
number  and  geometry  configuration  that  was  computed  previously  by  solving  the  steady 
Navier-Stokes  equations.  Static  pressure  distribution  evolution  along  the  centerline  in 
terms  of  a  pressure  coefficient  Cp  is  shown  in  Figure  45: 


It  can  be  seen  how  the  increasing  influx  of  kinetic  energy  is  converted  to  static 
pressure  as  the  solution  progresses  in  time.  After  sufficient  amount  of  time  the  flow 
transitions  from  separated  symmetric  to  separated  asymmetric  hence  the  drop  in  the 
overall  pressure  recovery.  The  red  lines  indicate  the  corresponding  steady  asymmetric 
solution  found  in  the  previous  chapter  for  the  same  geometry  and  Re  number.  Numbers  1 
through  6,  indicate  increasing  time. 
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2.  Flow  Visualization 

Velocity  and  streamline  contours  evolving  in  time  are  given  in  Fig.45  for  the 
examined  cases  below  the  lines  (a-a)  and  (b-b).  The  flow  starts  impulsively  and  ends  up 
to  a  steady  asymmetric  condition.  It  should  be  noted  that  before  the  asymmetric  flow 
pattern  is  reached  the  flow  temporarily  passes  through  a  symmetric  separated  condition 
and  stays  there  for  some  milliseconds  (time  steps)  before  numerical  truncation  errors 
perturb  the  flow  and  drive  it  to  the  final  more  stable  asymmetric  configuration. 
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Figure  47.  Typical  evolution  of  velocity  contours  in  time  for  a  case  below  the  flow 

regime  line  b-b. 
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For  cases  above  lines  a-a,  b-b  the  flow  starts  impulsively  and  after  sufficient  time 
it  does  not  approach  a  steady  condition.  On  the  contrary  the  flow  develops  initially  as 
separated  symmetric  but  after  some  time  it  has  been  found  to  break  down  shedding 
vortices  in  the  downstream  channel.  Eventually  these  unsteady  flow  structures  reach  the 
boundaries  where  they  no  longer  can  be  handled  computationally  because  the  boundaries 
downstream  and  upstream  are  invariant  with  time.  Thus  the  solution  accuracy  in  such 
cases  is  lost  a  few  time  steps  beyond  the  onset  of  flow  “break  down”.  Nevertheless  it  is 
indicated  that  for  cases  above  the  lines  a-a,  b-b  the  flow  seizes  to  be  two-dimensional.  It 
is  believed  that  the  flow  becomes  three-dimensional  and  unsteady  in  the  regime  above  the 
lines  given  in  Fig.29.  The  following  Figure  49  shows  velocity  contour  snap  shots  of 
successive  time  steps  of  such  a  case. 
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Figure  49.  Typical  evolution  of  velocity  contours  in  time  for  case  above  the  flow  regime 

line  b-b 
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VI.  CONCLUSIONS 


A  commercial  code  was  validated  against  the  analytical  solution  of  Jeffery-Hamel 
flows  through  two-dimensional  diffusers  with  less  than  1%  error.  The  code  was  found 
very  robust  and  user  friendly  and  was  used  with  a  great  degree  of  confidence  for  the 
subsequent  computations  of  diffuser  flow  regime  and  performance  maps.  A  large  number 
of  diffuser  geometries  (area  ratios  AR  1.15-5  and  lengths  over  inlet  height  L/W|  1-48) 
were  considered  to  compute  flow  regime  and  pressure  recovery  maps  for  incompressible 
laminar  flow  through  2-D  plane  diffusers  by  purely  computational  means  in  the  Re  range 
105-1048. 

A  threshold  combination  of  AR-L/Wi-Re  exists  beyond  that  the  flow  seizes  to  be 
steady  two-dimensional.  For  AR-L/Wi-Re  combinations  were  the  flow  remains  two- 
dimensional  both  symmetric  and  asymmetric  solutions  of  the  equations  are  possible.  The 
general  trend  is  that  increasing  the  divergence  angle  of  diffuser,  the  pattern  of  the  flow 
changes  from  attached  to  separated  symmetric  to  separated  asymmetric  and  finally  to  a 
non  two-dimensional  configuration. 

Performance  and  flow  regimes  are  a  function  of  diffuser  geometric  parameters 
and  Reynolds  number.  While  Reynolds  number  increases  (from  105  to  629)  maximum 
pressure  recovery  exists  for  decreasing  diffusion  half  angles  (from  ~1 2  degrees  to  -  3 
degrees)  and  always  in  the  two-dimensional  flow  regime. 

The  performance  maps  show  that  using  a  diffuser  to  connect  two  channels  of 
different  heights  can  be  advantageous  compared  to  having  a  sudden  expansion.  In 
addition  the  optimum  region  on  the  maps  shifts  to  lower  angles  and  higher  diffuser 
lengths  as  Reynolds  number  increases.  Limitations  on  the  complete  construction  of  the 
maps  are  imposed  in  this  study  by  the  fact  that  numerical  calculations  are  possible  only 
when  the  flow  is  two-dimensional  (i.e.  the  program  converges). 
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APPENDIX  A.  MATLAB  CODE  ON  THE  SOLUTION  OF  JEFFERY- 
HAMEL  DIFFUSER  FLOW  NON-LINEAR  ODE  . 


clear 

global  Re_r  a 

Re_r=1400;  %  reduced  Re  number  :  Re_r=(a*Re),  Re=?*V*r/? 

a=2*5*pi/360;  %  dif f fuser/Nozzle  half  angle (use  +  for  diff.  &  -  for  nozz.) 

% 

A=linspace (O,1,100O) ; 

x=-85. 0484725;  %  initial  curvature 

[T,Y]  =  ode45(0hamel_l,[A],[l  Ox]); 

figure  (1)  plot  of  the  nondimensional  vel.  profiles  in  polar  coord. 

plot(T,Y(:,l),-T,Y(:,l)) 

grid  on 


%  x  VALUES  (INITIAL  CURVATURES)  FOR  THE  SHOOTING 
%  TO  CONVERGE  QUICKLY  ANT  ACCURATELY: 


%- 

— Re*a- 

-##— 

— 

DIFFUSER - 

— M— 

- NOZZLE - 

% 

O 

M 

x= 

-2.005O8465 

M 

x=-2. 00508465 

% 

50 

M 

x= 

-3.5394176 

M 

x=-l. 121980 

% 

100 

M 

x= 

-5.8691811 

M 

x=-0. 640165 

% 

118 

M 

x= 

-6.8802251 

## 

x=-0. 526902 

% 

300 

M 

x= 

-18.9682230 

## 

x=-0. 0912577 

% 

500 

M 

x= 

-32.4654512 

M 

x=-0. 018555 

% 

800 

M 

x= 

-52.3003336 

M 

x=-0. 0025015 

% 

1400 

M 

x= 

-85.0484725 

M 

x=-0. 0000955 

function  dY=hamel_l  ( t,y) 
global  Re_r  a 
dy=zeros (3,1) ; 
dy(l)=y(2); 
dy(2)=y(3); 

dy(3)=  -  4. *a. A2. *y (2)  -2. *a. *Re_r. *  y(l).*y(2); 
return 
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APPENDIX  B.  SELF-SIMILAR  FLOWS  IN  CONICAL  DIFFUSERS. 


Assuming  purely  radial  flow  the  velocity  components  are: 

yF{6,<t>) 

ur= - 2 -  ,  U0=  0  ,  Ui=\) 

r 

Substituting  in  the  momentum  equations  we  get: 
r-momentum 


dur  |  U e_  du r_  |  fa,  uj+U^=  1  dp  |  W2u 

dr  r  SO  rsm.0  d</>  r  p  dr 

2u2F2  _  1  dp  v2  f  d^F  1  d2F  cos 0_dF^ 

r5  p  dr  r4yd92  sin 2  6  d(f>2  sin  0  dO  ^ 


(37) 

(38) 


9-momentum 

1  dp  2v  du 

0  = - — +  — r - ^ 

rp  dO  r2  dO 

dp  _  2 pv2  dF 
~dd~  r3  ~dO 

O- momentum 


(39) 

(40) 


1  dp  2v  dur 

- - — — | - L 

rpsmOd<j)  r2  sin2  0  d(j> 

dp  _  2 pv2  dF 
d(j)  r3  sin#  d(j) 


(41) 

(42) 


Take  the  derivative  of  (38)  with  respect  of  0,  (p  respectively  and  substitute  the 
pressure  gradient  from  (40)  and  (42)  to  get  (44)  and  (46)  respectively: 
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=>  (43) 


8 

2  v2F2~ 

1  d 

f  5(0 

1  8 
_i _ 

V 

( d2F  1  d2F  cos  0  dF  \ 

do 

r5 

p  dr 

1 80) 

1  do 

r4 

Kd02  sin 2  0  d<f> 2  sin 0  dO  J 

4FgF_6aF  |  a3F  |  1  d3F  cosfld2F  2cos Od^F  1  dF 

r  dO  SO  dO3  sin 2  0  d<f>2  dO  sin  0  dO2  sin3#  d(j)2  sin 2  0  dO 


(44) 


And, 


8 

2v2F2  " 

1  d 

r  dp  ^ 

8 

V 

(d2F  1  d2F  cos  OdF^ 

d<j> 

r5 

p  dr 

wJ 

_i _ 

d(j) 

r4 

|  | 

yd02  sin 1 0  d(j)2  sin 0  dO  J 

(45) 


4fgf_  6  gf,  g3F  !  1  d3F  |  cos#  d2F 
r  d(f>  sin  0  dO  dO2  d</>  sin 2  0  dtf>3  sin  0  dOd(j> 


(46) 


Equations  (44)  and  (46)  form  a  system  that  must  be  solved  to  get  the  velocity 
profiles  F(0,  </>) .  However  r  remains  a  parameter  in  the  equations  thus  there  is  no  unique 
solution  independent  of  r.  For  instance  a  symmetric  with  respect  to  r  bell  shaped  velocity 
profile  has  dF  I  d<f>=0 .  Equation  (46)  is  satisfied  and  equation  (44)  becomes  an  ODE 
where  r  is  a  parameter  again: 

A_f  d^_6  dF_  |  J3F  |  cos##2F  1  JF 
r  dO  dO  dO3  sin  0  dO2  sin2  0  dO 

Thus  for  point  source  flows  through  conical  diffuser  there  exists  no  similarity 
solution. 
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APPENDIX  C.  SUDDEN  EXPANSION  FLOW  COMPUTATIONS 


A.  INVISCID  STREAMLINE  SLOPE  COMPUTATIONS 


u=ones (100,100) ;  %  initial  guess 

vl=100;  %  inlet  velocity  (m/sec) 

rho=1.2;  %  density  (Kg/m3) 

pl=100e3;  %  inlet  pressure  (Pa) 

[m,n] =size (u) ; 

%  boundary  conditions 
%  upper  boundary 
u (  1  ,  1:50  )  =  50; 

u(l:81,  50  )  =  SO; 

u (  81  ,  51:100)  =  50; 

%  lower  boundary 
u  (  100,1:100  )  =  0; 

%  inlet 

u(  :  ,  1  )  =  flipudflinspace  (0,50,100)  1  )  ;| 

%  outlet 

u(81: 100,100)  =  flipud(linspace (0,50,20)  1  )  ; 

% 

for  niter=l:1200  %  §  of  iterations 
for  j=2:n-51 
for  i=2:m-19 

u(i,j)=.25*(  u(i,j+l)  +  u(i,j-l)  +  ... 

u(i+l,j)  +  u(i-l, j) ) ; 

end 

end 

for  j  =2 : n-1 
for  i=82:m-l 

u(i,j)=.25*(  u(i,j+l)  +  u(i,j-l)  +  ... 

u(i+l,j)  +  u(i-l, j ) ) ; 

end 

end 

end 
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figure  (1) 
u(l: 80,51: 100) =NaN; 

phsi=[u(l:99,: ) ;flipud(u(l: 100, : ) ) ];  %  stream  function  values 
contour (phsi, 60)  %  plot  of  stream  function  inside  the  contraction 
title ( 1  Streamlines  ') 
axis  off 

%  :  %  difference  of  stream  function  values 

parall_diffl=10O*norm(u(: ,l)-u(: ,2) ) 
parall_diff2=100*norm(u[81: 100,75) -u( 81: 10O, lOO) ) 
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B.  PSEUDO-INVISCID  PRESSURE  RECOVERY  FOR  A  SUDDEN 

EXPANSION 

For  incompressible  flow  the  sum  of  the  forces  on  the  control  volume  indicated  in 
Fig.49  reads: 


PA -PiA  +  TAMm  =(pA2U2)U2-(pAlUl)Ul  (47) 


Where  p,  U,  A  are  static  pressures,  velocities,  areas  at  stations  1 ,  2  respectively.  It 
is  reasonable  to  assume  that  the  pressure  at  station  1  is  applied  to  the  whole  face  area  at 
station  1  i.e.  Ax=  A2.  Then  using  continuity  U]AI  =U2A2  to  eliminate  velocities  and 
dividing  by  the  dynamic  head  at  station  1,  equation  (23)  becomes: 


c  _  Pi~ Pi 
p  0.5  pU, 


(48) 


Where  Cp  is  the  coefficient  of  pressure  recovery.  The  total  pressure  loss 
coefficient  CpT  for  the  same  expansion  flow  configuration  reads: 


p 

n  _1t\ 

pT 


p 

1  T  2 


0.5  pU, 


(49) 


Where  PT  denotes  total  pressures.  For  potential  flow  through  the  expansion  the 
ideal  recovery  is: 
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1 


2 


(50) 


Cp,‘-=1~7 


K  /  A 


Equations  24,  25,  26  are  plotted  together  in  Fig.50.  An  area  ratio  of  2  gives  the 
max  recovery  of  0.5. 


VS 


Figure  5 1 .  Pressure  recoveries  Cp,  CpT,  Cpp. 
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APPENDIX  D.  FLOW  REGIME  MAP  EXTRACTION  METHOD 
FROM  THE  NUMERICAL  DATA 


The  following  figures  show  the  actual  numerical  data  from  which  the 
interpolating  lines  of  Figure  29  are  drawn.  Power  laws  of  the  form  (51)  interpolate  the 
numerical  data  points  at  the  uppermost  comers  of  the  poly-lines  shown  in  figures  in  order 
to  construct  the  flow  regime  map: 


AR  =  eh 


Wx 


(51) 


For  each  flow  regime  line  the  exponents  a,  b  are  shown  in  Table  4  below: 


Table  4  Flow  regime  map  line  exponents. 


Re 

a 

b 

105 

0.1976 

1.4009 

210 

0.2479 

1.1688 

314 

0.2261 

0.8350 

420 

0.1878 

0.5956 

629 

0.1773 

0.4399 

1048 

0.1626 

0.2245 

Figure  52.  Flow  regime  map  data  for  Re=105. 
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Figure  55.  Flow  regime  map  data  for  Re=420 


Figure  56.  Flow  regime  map  data  for  Re=629. 
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Figure  58. 


Figure  59. 


APPENDIX  E.  DIFFUSER  PERFORMANCE  MAPS  FOR 
REYNOLDS  NUMBERS  OF  105,  210,  314,  420  AND  629 


Diffuser  performance  map  for  Re=629  with  the  respective  asymmetric  flow 

regime  line. 


Diffuser  performance  map  for  Re=629  with  the  respective  asymmetric  flow 


regime  line. 
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Figure  60. 


Figure  61. 


Diffuser  performance  map  for  Re=629  with  the  respective  asymmetric  flow 


regime  line. 


Diffuser  performance  map  for  Re=420  with  the  respective  asymmetric  flow 

regime  line. 
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Figure  62. 


Diffuser  performance  map  for  Re=629  with  the  respective  asymmetric  flow 

regime  line. 
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